Entanglement on mixed stabiliser states: Normal Forms and Reduction Procedures 
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The stabiliser formalism allows the efficient description of a sizeable class of pure as well as mixed quan- 
tum states of n-qubit systems. That same formalism has important applications in the field of quantum error 
correcting codes, where mixed stabiliser states correspond to projectors on subspaces associated with stabiliser 
codes. 

In this paper, we derive efficient reduction procedures to obtain various useful normal forms for stabiliser 
states. We explicitly prove that these procedures will always converge to the correct result and that these proce- 
dures are efficient in that they only require a polynomial number of operations on the generators of the stabilisers. 

On one hand, we obtain two single-party normal forms. The first, the row-reduced echelon form, is obtained 
using only permutations and multiplications of generators. This form is useful to calculate partial traces of 
stabiliser states. The second is the fully reduced form, where the reduction procedure invokes single-qubit 
operations and CNOT operations as well. This normal form allows for the efficient calculation of the overlap 
between two stabiliser states, as well as of the Uhlmann fidelity between them, and their Bures distance. 

On the other hand, we also find a reduction procedure of bipartite stabiliser states, where the operations 
involved are restricted to be local ones. The two-party normal form thus obtained lies bare a very simple bipartite 
entanglement structure of stabiliser states. To wit, we prove that every bipartite mixed stabiliser state is locally 
equivalent to a direct product of a number of maximally entangled states and, potentially, a separable state. 
As a consequence, using this normal form we can efficiently calculate every reasonable bipartite entanglement 
measure of mixed stabiliser states. 

PACS numbers: 03.67.-a,03.65.Yz 



I. INTRODUCTION 

The exploration of the properties of quantum entanglement 
is one of the main branches of quantum information theory 
dill IS- While a reasonably detailed understanding of two- 
qubit entanglement has been achieved, the entanglement prop- 
erties of higher-dimensional or multi-particle systems remain 
largely unexplored, with only isolated results |4]. This is 
largely due to the complexity involved in these investigations, 
which in turn originates from the tensor product structure of 
the multi-particle state space. This structure leads to an expo- 
nential growth in the number of parameters that are required 
for the description of the state. The same problem arises when 
one attempts to consider the time-evolution of a many-body 
quantum system or, say, of a quantum computation. Gener- 
ally, a significant part of the Hilbert space and consequently 
an exponential number of parameters are required to describe 
the quantum system at all times. 

One way of approaching this situation is to impose con- 
straints on the set of states and/or the set of operations that 
one is interested in without curtailing the variety of possible 
qualitative entanglement structures too much. In this context 
an interesting class of states that arises is that of stabiliser 
states I SlsHlal^llOll w hich, via the concept of graph states 
O El ElTIl ElEa Q3 El have some connection with 
graph theory. The feature of these states that allows for a more 
detailed study of their entanglement properties is the fact that 
an n-particle stabiliser state is determined as the joint unique 
eigenvector with eigenvalue +1 of a set of only n tensor prod- 
ucts of Pauli operators. This results in a very compact de- 
scription of the quantum state requiring only of order 0(n 2 ) 
parameters and therefore provides hope that a more detailed 



understanding of their entanglement structure can be obtained. 
Despite this simplification, stabiliser states still exhibit multi- 
particle entanglement and permit, for example, the violation 
of local realism Hal . 

The stabiliser formalism not only allows for the efficient de- 
scription of a certain type of quantum states, but also permits 
the efficient simulation of a restricted, but nevertheless inter- 
esting, class of time evolutions, namely those that respect the 
tensor product structure of the Pauli operators 1 19, 
Again, these simulations can be performed in polynomial time 
in the number of qubits, in stark contrast to the simulation of 
a general time evolution of an rt-qubit system, which requires 
an amount of resources that is exponential in n. 

The stabiliser formalism uniquely specifies the quantum 
state of an n-qubit system employing only polynomial re- 
sources. This alone, however, is not sufficient. It is also im- 
portant to be able to derive all relevant physical quantities, es- 
pecially those relating to entanglement, directly from the sta- 
biliser formalism. Indeed, having first to deduce the state ex- 
plicitly and then computing the property from the state would 
generally involve an undesirable exponential overhead in re- 
sources. While one can expect a direct approach to be pos- 
sible in principle, it is evident that detailed and explicit pre- 
sentations of algorithms to achieve these tasks in a systematic 
way and whose convergence is proven are of significant inter- 
est. In the context of entanglement properties some effort has 
recently been expended in this direction in 12311 - where, em- 
ploying sophisticated tools from group theory, the existence 
of a useful entanglement measure for multi-particle stabiliser 
states was demonstrated. 

The present work progresses further in a similar direction. 
Employing elementary tools we present a number of nor- 
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mal forms for pure and mixed stabiliser states, together with 
explicit and detailed descriptions of algorithms (including 
proofs of convergence) that allow the reduction to these nor- 
mal forms. In turn, these normal forms then permit us to com- 
pute any entanglement measure, overlaps between stabiliser 
states and various other quantities. Detailed descriptions of 
the algorithms are provided that should make it straightfor- 
ward to implement these algorithms in any programming en- 
vironment and we are able to provide a (/9-tested) suite of Mat- 
Lab programs on request. 

This suite can then form the basis for more detailed studies 
and further applications of the stabiliser formalism to a whole 
range of physical questions (see also Il20l0 . This will be re- 
ported on in a forthcoming publication. 

The stabiliser formalism also plays a central role in the field 
of quantum error correcting codes. Mixed stabiliser states (de- 
fined in Section II) are in one-to-one correspondence with pro- 
jectors on subspaces associated with stabiliser codes |5[- Al- 
though our normal forms and reduction procedures have been 
designed with applications to entanglement theory in mind, 
they might have a bearing on stabiliser codes as well. 

The present paper is organized as follows: In Section II the 
basic notations and conventions we use are introduced while 
Section III describes the elementary operations that will form 
the basis of all reduction procedures. 

The single-party normal forms are the topic of Sections IV, 
V and VI. Section IV deals with the so-called row-reduced 
echelon form (algorithm RREF), the reduction to which is 
based on row operations only. It allows to check independence 
of any (putative) set of generators and to calculate partial 
traces (algorithm PTRACE). In Section V we describe the full 
reduction procedure (algorithm CNF1) to single-party normal 
form, using row and column (qubit) operations. In Section VI 
we present an algorithm (algorithm OVERLAP) that is based 
on the full reduction and allows to calculate overlaps between 
stabiliser states, Uhlmann fidelity, and Bures distance. 

In Section VII we turn to the bipartite case, where we prove 
that the bipartite entanglement structure of stabiliser states is 
remarkably simple. To wit, we show that mixed bipartite sta- 
biliser states are locally equivalent to a tensor product of a cer- 
tain number of maximally entangled 2-qubit states and, poten- 
tially, a fully sparable mixed state. We present an algorithm 
(algorithm CNFP) to obtain the number of these maximally 
entangled pairs, allowing for the calculation of any reasonable 
bipartite entanglement measure. 

We conclude the description of our findings in Section VIII. 

II. NOTATIONS AND CONVENTIONS 

A stabiliser operator on N qubits is a tensor product of 
operators taken from the set of Pauli operators 

*=-(!i)."=-(!7).*=-(i-°i). a > 

and the identity 1. An example for N — 3 would be the 
operator g — X ® 1 ® Z. A set G — {gi, . . . , gx } of K 
mutually commuting stabiliser operators that are independent, 



i.e. Yli=i 9i* = ^ exactly if all are even, are called the 
generator set for the stabiliser group S. This stabiliser group 
S then consists of all distinct products of operators from the 
generator set. 

For K — N a generator set G uniquely determines a single 
state \ip) that satisfies gk\ip) = \'4 ) ) f° r all fc = l,...,N. Any 
state for which such a generator set exists is called stabiliser 
state. Such a state has trivially the property that gk\ip)(if>\ = 
\i>)(ip\ for all k so that 

MM = i5> (2) 

This formula depends on the complete set of stabilisers, and is, 
therefore, not very practical. The following formula expresses 
the state in terms of a generator set 1 8] 

N 

fe=i 

The procedure presented in Section ffVl yields as a side result 
an elementary proof of this statement. 

Considering two parties A and B, the reduced density ma- 
trix of the stabiliser state can be computed as 

p A = Tr B p=^ w J2 Tr B9- (4) 

This obviously means that only operators g contribute that 
have identity operators acting on all qubits belonging to B. 
Needless to say, computing pa from p directly is hopelessly 
inefficient; as there are 2 N different g, this task requires an 
exponential number of operations in general. 

It turns out, however, that there is a class of mixed states 
that can also be characterised employing stabilisers. We will 
call these mixed stabiliser states, and they contain the pure 
stabiliser states as a subset. The important feature of this class 
is that the reduced density matrix of a mixed stabiliser state is 
again a mixed stabiliser state. Furthermore, as we will show 
below in Section lTV CI the stabiliser group of a reduction can 
be efficiently calculated directly from the original stabiliser 
group, without calculating the state and its reduction explic- 
itly. 

To characterise mixed stabiliser states, one simply consid- 
ers sets G that are linearly dependent. As a consequence, by 
multiplying stabilisers, one can achieve that some of them be- 
come identical to 1 and only K linearly independent ones re- 
main. Then the common eigenspace of these K operators will 
have a dimension larger than 1. As in eq. Q one immediately 
deduces that the density operator is again just the projector 
onto this eigenspace, rescaled to trace 1 : 

v = ^1 + ^ 

fc=i 



Given that V is a projector onto a subspace of dimension 
2 N ~ K , the entropy of p is simply N — K. In analogy with 
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matrix analysis terminology, we will call K the rank of the 
stabiliser group. Stabiliser groups with K = N will be called 
full-rank, and stabiliser groups with K < N rank-deficient. 

Note: In case of a rank-deficient stabiliser group one has 
to distinguish between stabiliser states and stabilised states. 
The stabiliser state is the one given by (0, and the stabiliser 
formalism allows to study its properties in an efficient way. 
On the other hand, there are many states that are stabilised 
by that same stabiliser group, but in general they are not sta- 
biliser states. Indeed, most of these stabilised states cannot 
be described as "the unique state stabilised by a full-rank sta- 
biliser group", and hence, the stabiliser formalism cannot be 
used to study their properties via that group. For example, 
any state is stabilised by the (singleton) stabiliser group {1}, 
but only the maximally mixed state 1 /2 is the stabiliser state 
corresponding to that group. 

For the following the aim will be to derive basic entangle- 
ment properties such as the entropy of entanglement or the 
logarithmic negativity for stabiliser states, pure or mixed, di- 
rectly from their generating set. To this end it will be useful 
to find a normal form for the generator set that reveals the 
relevant entanglement structure. 

We now introduce the concept of stabiliser array, which is 
quite simply a rectangular array of K rows and N columns, 
where the elements are Pauli matrices or the Identity ma- 
trix. Specifically, the element in the fc-th row and n-th col- 
umn of the stabiliser array corresponding to a generator set 
G = {g%, . . . , gx} on N qubits is the n-th tensor factor (cor- 
responding to qubit n) of the fc-th generator gk . For some ap- 
plications it will also be necessary to deal with the generator 
phase factors. While in general these phase factors can assume 
the values ±1 and ±i, for the purpose of describing stabiliser 
states only the values ±1 make sense (since states are Hermi- 
tian). We will store these phase factors in a if-dimensional 
vector s, where sj. is the phase factor of generator gk- 

The purpose of the various normal forms that will be pre- 
sented in this paper are to structure the set of stabiliser states 
into certain equivalence classes. They are similar in spirit to 
the normal forms that have been devised for matrices. For ex- 
ample, the row-reduced echelon form, which exhibits the rank 
of a matrix, has a counterpart for stabiliser arrays. Despite this 
similarity, the normal forms presented here are of an entirely 
different nature. The rank of a stabiliser array, which we will 
introduce in Section HVl is akin to the rank of a matrix in that 
it equals the number of independent generators in a generator 
set, but there the similarity stops. In linear algebra one de- 
fines both row and column rank of a matrix and one proves 
that these two ranks are actually equal. For stabiliser arrays, 
one cannot even give a meaningful definition of column rank. 

These differences ultimately boil down to the fact that a sta- 
biliser array is not really a matrix. The two foremost reasons 
are that its elements are not numbers but elements of the Pauli 
group, and second, that matrices represent linear operations 
in linear spaces, while stabiliser arrays represent sets (namely, 
sets of generators). As a consequence, while operations like 
matrix transpose, matrix multiplication, addition, and inverse 
make perfect sense for matrices, they are utterly meaningless 
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TABLE I: Multiplication table for Pauli operators; shown is a row .a co \. 



for stabiliser arrays. The allowed operations on stabiliser ar- 
rays are thus much more restricted than in the matrix case. 
For example, the only row operations that make sense for sta- 
biliser arrays are row interchange and elementwise row mul- 
tiplication (which is based on the Pauli group multiplication 
law). This will be discussed in more detail in the following 
Section. This explains the need for entirely new reduction 
procedures for stabiliser arrays. 



m. ELEMENTARY OPERATIONS 

In this Section, we describe the allowed elementary opera- 
tions that transform a stabiliser array and which we will use 
to reduce an array to its normal forms. As in the matrix case, 
these operations come in two kinds. The first kind are the row 
operations. It is important to realise that row operations will 
not alter the stabiliser state at all, but only alter the genera- 
tor set it is represented by. These are the row transposition, 
which interchanges (transposes) two rows in the stabiliser ar- 
ray, and the row multiplication, which multiplies one row with 
another one. The latter operation changes the generators of the 
stabiliser group, but not the group itself and hence not the sta- 
biliser state either. We will use the phrase "multiply row fc 
with row I" to mean "multiply rows k and I elementwise and 
set row I to the product obtained." The multiplication table for 
Pauli operators is shown in Table U 

The second kind of operations are the column operations, 
which may alter the state. The column operations we will use 
are a certain class of single-qubit operations, transposing two 
columns, and the CNOT operation between two qubits. As 
single-qubit operations we take those that act on one given 
column of the stabiliser operators by permuting the Pauli op- 
erators (in the given column) among themselves. These op- 
erations can be constructed from combinations of Hadamard 
gates (H) and n /4 gates (P) (see Table [n}. Note that odd 
permutations must involve a sign change in one of the Pauli 
operators in order to correspond to a unitary operation. The 
particular sign changes of Table UTI have been chosen to make 
the unitaries implementing the odd permutations involutory 
(apart from a global phase). That is, UU = exp(icf))l. Note 
also that the second and third permutation in the Table are 
each other's inverse. 

For the bipartite normal form described in Section IVTT1 we 
will need to divide the qubits into two parties and only al- 
low operations that are local to those parties. Transposing 
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in a recursive fashion. There are three cases: 



TABLE II: Truth table for the single-qubit operations employed by 
the CNF algorithm. Any permutation of the set of Pauli operators 
can be achieved. 
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TABLE III: Truth table for the CNOT gate employed by the CNF 
algorithm. C and T refer to control and target qubit, respectively. 
The primed columns give the values after the operation. 



two columns in the bipartite case is only allowed when both 
columns (qubits) belong to the same party. Otherwise this 
would be a non-local operation, which would very likely af- 
fect the amount of entanglement. 

The CNOT gate between two qubits, one being the control 
qubit and one the target qubit, operates on the two correspond- 
ing columns of the stabiliser array. In the case of a bipartite 
CNF, we must again ensure that both qubits belong to the same 
party. The truth table for the CNOT is given in table [HI] Note 
that the column pertaining to the control qubit is modified too; 
this is a peculiarity of the description of states by stabilisers. 
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The symbols ':' and '. . .' denote a number of repeated rows 
and columns. This number may be zero. The symbol RREF' 
denotes a (possibly empty) sub-array that is also in RREF 
form. The symbol * denotes either a Pauli operator or an iden- 
tity 1 . Furthermore, a, o\ and 02 are Pauli operators, and <j\ 
and <7 2 anticommute. We will refer to the operators in these 
positions as column leaders of their column, and row leaders 
of their row. 

The RREF algorithm works by applying a sequence of el- 
ementary row operations to the stabiliser array. At every step 
of the algorithm it is determined which elementary operation 
to apply based on the values contained in a certain contiguous 
subarray of the full array. At every step this subarray, which 
we will call the active region, either stays the same or de- 
creases in size. The algorithm terminates when the size of the 
active region has decreased to zero. Note that the elementary 
operations operate on the full stabiliser array and not just on 
the active region. 

Let K, N be the number of rows (generators) and columns 
(qubits) of the stabiliser array, respectively. The variable Kjj 
contains the index of the first row in the active region, and Nl 
the index of its first column. The active region thus consists of 
the array elements for Kjj < i < K and Nl < j < N. 
Initially, the active region comprises the full stabiliser array, 
hence Ky = 1 and Nl = 1. 

In this and subsequent sections, the phase factors exp(i^fc) 
of the various generators will not be mentioned explicitly. 
They are best maintained under the form of a single additional 
column in the stabiliser array, which is modified by row per- 
mutations and the elementary operations of tablesHlllTlandlllll 



A. Algorithm RREF 



IV. ROW-REDUCED ECHELON FORM 

While the Clifford Normal Form (CNF) of a stabiliser array 
will be obtained below via application of both elementary row 
and column operations, it is possible to obtain a normal form 
using elementary row operations only. Due to its similarity 
to the matrix case, we will call this normal form the Row- 
Reduced Echelon Form (RREF). The benefits of the RREF are 
that it is very easy to obtain, the stabiliser state represented by 
the stabiliser array is not changed, and it is applicable to states 
on any number of parties. Furthermore, as we shall see below, 
it is an efficient way to eliminate linearly dependent rows from 
the stabiliser array. 

The general structure of the RREF is most easily described 



1. Count the number of different Pauli operators (X, Y 
and Z) in the first column (Nl) of the active region, i.e. 
restricting attention only to rows Ku up to K. 

2. Three cases can be considered: 

(a) There are no Pauli operators in column Nl- 

i. Increase Nl by 1 . 

(b) There is only 1 kind of Pauli operator. 

Let k be the first row in the active region where 
column Nl contains a Pauli operator. 

i. Make row k the top row of the active region 
by transposing, if necessary, row k with row 
K v . 
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Initial stabiliser array: 




Depending on the content of row 3, do the following: 
1 : Do nothing. 
<Ti: Multiply row 1 with row 3. 
(72 : Multiply row 2 with row 3. 

173: Multiply row 1 with row 3, and then row 2 with row 3. 



having rik leading 1 operators, and having Ok as row lead- 
ers. Let the generators be sorted according to nk in ascending 
order. When multiplying two rows that satisfy 712 > n\, the 
number of leading 1 operators in the product is n\ and the row 
leader is either o\ (if 77,2 > n{) or o\o~2 (if "2 = "-1), which is 
different from either o\ and 02- In both cases this shows that 
the product cannot occur as another generator in the RREF. 
This proves that it is not possible to write one RREF genera- 
tor as a product of other RREF generators. □ 



C. Partial Trace of a Stabiliser State 



TABLE IV: Required operations to eliminate any Pauli operator from 
row 3 of the stabiliser array shown above. The operators ai, <J2, and 
03 are a permutation of X, Y and Z. 



ii. Multiply row Kjj with all other rows in the 
active region that have the same Pauli in col- 
umn Nl- 

iii. Increase Kjj and Nl by 1. 

(c) There are at least 2 different kinds of Pauli opera- 
tors. 

Let k\ be the first row in the active region where 
column Nl contains a Pauli operator, and fc 2 be 
the first row in the active region where column Nl 
contains a different Pauli operator. 

i. Make row k\ the top row of the active region 
by transposing, if necessary, row k\ with row 
K v . 

ii. Make row k% the second row of the active re- 
gion by transposing, if necessary, row fc 2 with 
row Kjj + 1, 

iii. Multiply every other row in the active region 
with either row Kjj, row Kjj + 1, both rows 
or none, depending on the element in column 
Nl (see TablelFvl 

3. If the active region still has non-zero size {Nl < N and 
Kjj < K), continue with step 1, else terminate. 



B. Checking independence of a set of generators 

The easiest way to check independence of a set of genera- 
tors is to compute the RREF of the stabiliser array. This fact is 
one other property the stabiliser RREF and the matrix RREF 
have in common. Dependencies between generators will show 
up as RREF rows containing only 1 operators. Removing 
these all- 1 rows leaves an independent set of generators. 

Proof. From the form of the RREF one observes that there 
cannot be more than two rows with the same number of lead- 
ing 1 operators, and if there are two such, they have a different 
row leader. Consider a subset of generators gt in the RREF, 



A useful and important operation is the partial trace. The 
RREF algorithm is the central part in the following efficient 
partial trace algorithm: 

Algorithm PTRACE 

1. By column permutations bring the columns of the 
qubits to be traced out in first position. 

2. Bring those columns to RREF. 

3. Remove the rows containing the column leader(s). 

4. Finally, remove those columns themselves. 

Proof. To prove that this algorithm indeed calculates the 
partial trace, consider again the three cases for the RREF: 
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We have to show that the state described by RREF', say p', is 
the state obtained from the original stabiliser state p by tracing 
out the qubit pertaining to column 1 . Denote the sequences of 
* operators by g, g\ and 172, respectively. 

Using eq. l|5}, it is easy to see that, in the first case, 



1 

P= 2 



1 P 



and tracing out qubit 1 yields 
In the second case, 

1 <g) 1 + a ® g 



P 



(1 <8> p) 



-(1 ®p' + <r®gp'), 



and again, as Pauli operators have trace 0, 
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In the third and final case, 



completion of the algorithm, is 



P = 



1 + <J\ <g> gi l®t+<T2®92 



resulting yet again in 



T ri p = p'. 



(1 <8> 2/9') 



+ (Ti(T2 ® gi92p')i 



□ 



/X 1 

1 X 



1 1 



1 1 ... X 



1\ 

1 



1 ... 1 



1/ 



Here we have left open the possibility that the rows of the 
initial stabiliser array might not be independent. 



A. Algorithm CNF1 



V. SINGLE-PARTY NORMAL FORM 

The CNF algorithm works by applying a sequence of el- 
ementary operations to the stabiliser array. At every step of 
the algorithm it is determined which elementary operation to 
apply based on the values contained in a certain contiguous 
subarray of the full array. At every step this subarray, which 
we will call the active region, either stays the same or de- 
creases in size. The algorithm terminates when the size of the 
active region has decreased to zero. Note that the elementary 
operations operate on the full stabiliser array and not just on 
the active region. 

Let K, N be the number of rows (generators) and columns 
(qubits) of the stabiliser array, respectively. The variables 
Kjj and Kl contain the indices of the first (uppermost) and 
last (lowermost) row in the active region, and N L and N R 
the indices of its first (leftmost) and last (rightmost) column. 
The active region thus consists of the array elements for 
Kjj < i < Kl and Nl < j < Nr. Initially, the active region 
comprises the full stabiliser array, hence Ku = 1, Kl = K, 
N L — 1 and N R — N. We will prove below that after ev- 
ery iteration of the algorithm the stabiliser array has the block 
structure 



fx 1 

1 X 



1 1 



X 



1 1 



V i i 
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1 ... 1 
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The block containing the asterisks is the active region and has 
not yet been brought to normal form. The columns on the left 
of the active region correspond to qubits that are in an eigen- 
state of the X operator, the columns on its right correspond to 
qubits that are in a totally mixed state. The final form, after 



1. Count the number of different Pauli operators (X, Y 
and Z) in the first column (Nl) of the active region, i.e. 
restricting attention only to rows Kjj up to Kl- 

2. Three cases can be considered: 

(a) There are no Pauli operators in column Nl- 

i. If necessary, transpose column N L with col- 
umn Nr. 

ii. Decrease Nr by 1. 

(b) There is only 1 kind of Pauli operator. 

Let k be the first row in the active region where 
column Nl contains a Pauli operator. 

i. Make row k the top row of the active region 
by transposing, if necessary, row k with row 
Ku. 

ii. Apply whatever single-qubit operation on 
column N L that brings that Pauli operator to 
an X. 

iii. Multiply row Ku with all other rows in the 
active region that have an X in column Nl- 

iv. Consider the elements of the first row of the 
active region (row Ku). To each of the 
columns beyond the first one that contains in 
the first row a Pauli different from X, apply a 
single-qubit operation to turn it into an X. 

v. To each of these columns, which now have an 
X in the first row, successively apply a CNOT 
operation with control column Nl- 

vi. Increase Ku and Nl by 1 . 

(c) There are at least 2 different kinds of Pauli opera- 
tors. 

Let ki be the first row in the active region where 
column Nl contains a Pauli operator, and ki be 
the first row in the active region where column Nl 
contains a different Pauli operator. 

i. Make row k\ the top row of the active region 
by transposing, if necessary, row k\ with row 
K v . 
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ii. Make row &2 the second row of the active re- 
gion by transposing, if necessary, row fc 2 with 
row Ku + 1. 

iii. Bring the element on row Kjj to an X and 
the element on row Kjj + 1 to a Z by apply- 
ing, if necessary, a single-qubit operation on 
column Nl ■ 

iv. Consider the first two rows of the active re- 
gion (rows Ku and Ku + 1). Find the first 
column beyond column Nl, say column I, 
that contains an anticommuting pair on those 
rows (i.e. two non-identical Pauli operators). 

v. Bring the anticommuting pair to an (X, Y) 
pair by applying, if necessary, a single-qubit 
operation to that column. 

vi. Apply a CNOT operation to that column, 
with column Nl as control. 

vii. The extent of the active region is not changed 
in this case. 

3. If the active region still has non-zero size (Nl < Nr 
and Ku < Kl), continue with step 1, else terminate. 



B. Proof of correctness of algorithm CNF1 

We will now show that algorithm CNF1 indeed brings any 
stabiliser array into its normal form. We consider the three 
cases (a), (b) and (c) in succession. 



1. Case (a) 

This case corresponds to column Nl containing 1 only and 
therefore belongs to the block right of the active region. Step 
(a.i) does just that and step (a.ii) subsequently excludes this 
column from the active region. 

2. Case (b) 

This case corresponds to a column containing 1 operators 
and Pauli operators of just one kind. Step (b.i) brings the first 
of these Pauli operators to the top row, with the ultimate goal 
of excluding this row from the active region. Step (b.ii) applies 
a single-qubit rotation to bring the Pauli operators in standard 
form, which in this case is an X operator. 

In step (b.iii) the column is then "cleaned up". Through 
multiplying the top row Ku with other rows containing an X 
in column Nl, we obtain a stabiliser array that is still describ- 
ing the same state but contains only one X in column Nl- So 
this column is already in standard form. 

However, the top row is not in standard form yet. Step (b.iv) 
applies an appropriate single-qubit operation to every column 
in the active region, except for the first one, so that the first row 
contains either 1 or X operators. Step (b.v) then performs a 
"row cleanup", by applying CNOTs to the columns starting 



with an X, the first column being the control column. The 
target X operators are thereby turned into 1, leaving us with 
a first row of the form (X, 1 , . . . , 1 ) . 

It is not a priori clear, however, that step (b.v) is not undo- 
ing the cleanup of column Nl by step (b.iii). Nevertheless, 
inspection of the CNOT truth table reveals that the 1 opera- 
tors in column Nl can either be turned into a Z or remain 1, 
by any number of CNOTs. Although a Z operator can actu- 
ally occur during the execution of step (b.iii), in the end all 
operators will be turned back into 1. This must be so because 
the top row of the active region is turned into (X, !,...,!), 
which does not commute with a row starting with a Z. So the 
assumption of commutativity of the generators ensures that 
step (b.iii) is not undone by step (b.v). 

Finally, we note that both the first row and the first column 
are now in standard form and can be removed from the active 
region, which is done in step (b.vi). The top left block in the 
normal form array hereby receives one further X operator. 



3. Case (c) 

The most difficult case to investigate is case (c), because 
here it is not a priori clear that any progress is made within an 
iteration. Indeed, the extent of the active region is not changed 
and it is not clear that further iterations will eventually escape 
from case (c), resulting potentially in an infinite loop. 

However, every execution of case (c) does result in mea- 
surable progress. As can be seen from the truth table of the 
CNOT operation, the end result of the CNOT in step (c.vi) is 





N l 


I 




I' 


Ku 


X 


X 


X 


1 


Ku + 1 


Z 


Y 


t 


Y 



Hence, a 1 is introduced in row Ku where there originally 
was none. Furthermore, no further algorithmic step in case (c) 
ever touches this element again, by the following reasoning. 

• The only operations that do change the top row Ku are 
the transposition in step (c.i), and the CNOT in step 
(c.vi). 

• Step (c.i) is executed at most once before the algorithm 
breaks out of the (c) case, namely at the very beginning. 
This is because the X brought in the top left position is 
not changed by the CNOT. 

• The first CNOT that operates on target column I intro- 
duces the 1 there. In further iterations, the CNOT will 
not operate on column I a second time, because step 
(civ) sets the target column to a column containing an 
anticommuting pair in the top two rows, and the 1 cre- 
ated in the top of column I does not form part of an 
anticommuting pair. 

It is now easy to see why the algorithm must eventually break 
out of the (c) case. Every iteration through this case increases 
the number of 1 operators in the top row by 1, but there are 
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only a limited number of places (columns) available to do this. 
Hence the number of successive iterations through case (c) 
must be limited too. 

The algorithm breaks out of the loop through the (c) case 
when there are no further anticommuting pairs in column Nl . 
As a consequence, the algorithm will then either execute case 
(a) or case (b), thereby again reducing the extent of the active 
region. 

C. Alternative proof of projection formulas and {5j 

In this subsection we present a new proof of the equivalence 
of the expressions (0 and for a pure stabiliser state, and of 
(0} and (0 for a mixed stabiliser state. 

By the proof of the CNF1 algorithm, a state described by a 
certain stabiliser array is unitarily equivalent to the state de- 
scribed by the normal form of that array. Let the initial sta- 
biliser group S be given by a stabiliser array. Let S' be the 
stabiliser group described by the normal form of that array. 
The K generators of S 1 are of the form 

g' k = 1 ® . . . ®X® . . . ® 1, 

with the X operator in the fc-th tensor factor. The stabiliser 
state corresponding to the normal form is therefore 

P ' = x h ®...®x iK ®i® N - K 

ii,...,i K e{0,l} 

= ((1 + X)/2)® K ® (1/2)® N ~ K 
= 1 TT 1+ gfc 

2 N-K 11 2 ' 
fe=l 

Let U be the unitary corresponding to the sequence of elemen- 
tary operations that brought the stabiliser array to its normal 
form. To wit, S consists of the elements g = Ug'W, g' € S', 
and can be generated by the generators gk := Ug k W . Then 
the stabiliser state corresponding to S is given by 

P = ^nY,9 

ges 

= ^ E 

g'es' 
= Up'U^ 

2 N-K LL U 2 

fc=l 

1 TT 1 +9k 

2 n-k TT 2 ■ 

k=l 

□ 



VI. FIDELITY BETWEEN STABILISER STATES 

The topic of this section is an algorithm to calculate the 
overlap F = Tr[pxp2] between two mixed stabiliser states p\ 



and p2 directly from their K± x N and K% x N stabiliser 
arrays A\ and Ay,. 

While the overlap between two states is certainly an inter- 
esting quantity, the Bures distance 

D( Pl ,p 2 ) := 2 y /l-F u (p 1 ,p 2 ), 
where F u is the Uhlmann fidelity 

F u (fil,p2) ■■= Tr[Y y/PlP2y/pi\, 

is a much more desirable quantity, as it is an actual distance 
measure and has a much nicer interpretation. It is well-known 
that for pure states the Uhlmann fidelity between two states is 
just the square root of their overlap, while for general mixed 
states there is no such relation. It will turn out that with just 
a minor modification the algorithm is also able to calculate 
the Uhlmann fidelity. This allows us to calculate the overlap, 
the Uhlmann fidelity and Bures distance for stabiliser states in 
one go. 

For the calculation of the overlap (and fidelity) it is impera- 
tive to take the generator phases into account. We will use the 
vectors Si and S2 for that purpose. The elementary row op- 
erations of multiplication and permutation of stabiliser rows 
are understood to treat the phase vector as an additional col- 
umn of the stabiliser array. Furthermore, row multiplication, 
single-qubit rotation and CNOT operation have to multiply the 
appropriate generator phase with the phase factor mentioned 
in their truth tables. 



A. Algorithm OVERLAP 

1. Construct the (K\ + K2) x N composite array A and 
the composite vector S of generator phases: 

'-(£)■ 

124] 

2. By applying the CNF1 algorithm to the composite array 
A (and its vector S of generator phases), with initial 
active region set to the full A\ part (excluding A2 !), the 
A\ part is brought to CNF form, while automatically 
applying the same sequence of column operations to the 
A 2 part. Let Ri be the number of X operators in this 
CNF. 

3. Set the active region to all the rows of the A 2 part and all 
the columns for which A\ contains X operators. That 

is, Kjj = Ki + 1, K L = Ki + K 2 , N L = 1 and 
N R = R. Set T = 1. 

4. Count the number of different Pauli operators (X, Y 
and Z) in the first column (Nl) of the active region, i.e. 
restricting attention only to rows Kjj up to Kl- 

5. Three cases can be considered: 
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(a) There are no Pauli operators in column Nl- 

i. Do nothing. 

(b) There is only 1 kind of Pauli operator. 

Let k be the first row in the active region where 
column Nl contains a Pauli operator. 

i. Make row k the top row of the active region 
by transposing, if necessary, row k with row 
K v . 

ii. Multiply row Kjj with all other rows in the 
active region that have a Pauli operator (nec- 
essarily equal to the one on row Kjj) in col- 
umn Nl- 

iii. Let P be the element in column Nl on row 
Kfj (the column leader). 

• If P is not an X, divide T by 2 and in- 
crease Kjj by 1 . 

• If P is an X, multiply row Nl (that is, 
the row containing an X in column Nl 
of subarray Af) to row Kjj. 

(c) There are at least 2 different kinds of Pauli opera- 
tors. 

Let k\ be the first row in the active region where 
column Nl contains a Pauli operator, and &2 be 
the first row in the active region where column Nl 
contains a different Pauli operator. 

i. Make row k\ the top row of the active region 
by transposing, if necessary, row ki with row 
K v . 

ii. Make row k 2 the second row of the active re- 
gion by transposing, if necessary, row k 2 with 

row Kjj + 1, 

iii. Multiply every other row in the active region 
with either row Kjj, row Kjj + 1, both rows 
or none, depending on the element in column 
Nl (according to Table llVl i. 

iv. Let Pi be the element in column Nl on row 
Kjj, and P2 the one on row Kjj + 1 (the col- 
umn leaders). Turn P2 into an X, as follows: 

• If neither Pi nor P2 is an X, multiply 
row Kjj to row Kjj + 1, effectively turn- 
ing P2 into an X. 

• If Pi is an X, transpose row Kjj with 
Ku + 1. 

v. Multiply row Nl (that is, the row containing 
an X in column Nl of subarray Af) to row 
Ku + 1. 

vi. Divide T by 2, and increase Kjj by 1. 

6. Increase Nl by 1. 

7. If the active region is still non-empty {Nl < Nr and 
Kjj < Kjf), continue with step 4. 

8. (End Game) Here we calculate a correction factor C for 
the overlap and the fidelity. Set C — 1 as default value. 
If Kjj < Kl do the following. 



(a) Case Nl < N: Consider the bottom right block 
of stabiliser array A2 consisting of rows Kjj to 
Kl and columns Nr + 1 to N. Calculate the rank 
i?2 of that block using, e.g. the PvREF algorithm. 
From i? 2 calculate the correction factor as C = 

(b) Case Nl > N: Let tk be the generator phase of 
row k. If at least one of the tk for Kjj < k < Kl 
is -1, setT = 0. 

Terminate with return values F = CT/2 N - Kl+N ~ K2 
for the overlap and F u = C y/T /2 N -^+ N - K ^ for the 
Uhlmann fidelity. 



B. Proof of correctness of algorithm OVERLAP 

The overlap F — Tr^i^] can be calculated iteratively by 
performing the trace as a succession of partial traces over sin- 
gle qubits: F = Tr[Tri [,01/02]], where Tri denotes the par- 
tial trace over the first qubit. What we need to show is that 
one iteration of steps 4-7 indeed performs this single-qubit 
partial trace. It will be convenient to express the overlap in 
terms of the projectors Vi and V2, with p\ = Vi/2 N ~ Kl and 
92 = V 2 /2 N - I<2 . Then 



1 



2N-K 1 +N-K 2 



Keeping in mind that we also want to calculate F u , we will 
proceed by first calculating V\PiV\. The overlap is just the 
trace of this quantity, which is the same as TrfPiT^] by virtue 
of Vi being a projector. 

Step 2 of the algorithm applies the same sequence of uni- 
taries to both states, hence the overlap between them does not 
change (and neither does the Uhlmann fidelity). Let Vi thus 
be specified by a CNF stabiliser array, containing Ri < N 
X -operators: 



Vi = 



1 + SiX 
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If all generators in array A\ are independent, we obviously 
have Ri — Ki . In the above expression, Si is the generator 
phase of the i-th generator of Vi . Likewise, we will denote 
by ti the generator phase of the z-th generator of V 2 - Further- 
more, let V[ be the stabiliser projector of the array obtained 
by deleting row 1 and column 1 from A\. 

In the following, we will calculate V1V2V1 and show that 
it is equal to a certain scalar value T times a tensor product 
of rank-1 projectors and identity operators. We will proceed 
in an iterative fashion, by showing that ViP%P\ decomposes 
as a scalar T\ times either a rank-1 projector or an identity 
tensored with a product V[V 2 Vi of projectors over qubits 2 to 
N. To calculate T, we start off with the initial value T = 1 
and update T by multiplying it with the value of Ti found at 
each iteration. 

We will assume first that A\ and A2 have more than 1 col- 
umn. The case that they only have 1 column, which is what 
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can happen in the final iteration of the algorithm, will be con- 
sidered in subsection 4. We will also assume that the first ten- 
sor factor of V\ is a rank-1 projector (i.e. R\ > 0). The case 
that Vi equals the identity (which will again typically happen 
at the end of the iterations) will also be covered in subsection 
4. 

Let us now take on the main case, where there are at least 
two tensor factors to consider, and the first factor of V\ is 
t+ * iX . Thus we can write V\ = t+ ° iX ®V[. As in algo- 
rithm CNF1 there are three cases to consider, depending on 
the number of different Pauli operators contained in the first 
column of the second array A 2 . We will investigate these three 
possibilities in succession. 

It is useful to note that 



1 +sX 1 
z cr — 



sX 



(l+sX)/2, £7 = 1 

s(l+sX)/2, a = X 
0, a = Y,Z 



We can therefore distinguish two cases. If a is not an X, we 
find 

VxV 2 Vx = \ t + * lX ®PiV 2 Vi. 

This corresponds to a value of T\ = 1/2. This is implemented 
in step (b.iii, first case) by dividing the running T by 2, and 
deleting row 1 and column 1 from A 2 . 
If, on the other hand, a = X, we have 

v 1 v 2 v 1 = \±±^L ®v[v' 2 v[ 

Siti 1 + S\X . 
+ ^ ^—®V[g'V 2 V' x 

<g> V[ V 2 V[ 



2 

1 + Sl X 



1. Case (a) 



where V 2 = ^ +s ^ tl9 ~p' 2 is a projector. This corresponds to 
a scalar value of T\ = 1. This is accomplished in step (b.iii, 
second case) by multiplying row 1 of A\ to row 1 of A 2 , and 
deleting column 1 of A 2 (leaving row 1). 



If the first column of A 2 contains no Pauli operators, this 
corresponds to V 2 being of the form 

V 2 = t®V 2 , 

where V 2 is the stabiliser projector of the array obtained by 
deleting column 1 from A 2 . Hence 

ViV 2 Vi = ( 1 + * lX v[) (i ® v 2 ) ( 1 + ° lX ® n) 



1 + siX 



<V[V 2 V[ 



This is indeed of the form claimed above, with scalar value 
T\ = 1. Hence, nothing needs to be done in this iteration 
except for deleting column 1 . 



2. Case (b) 

Steps (b.i) and (b.ii) bring column 1 of V 2 to RREF form. 
In this case, column 1 will contain a single Pauli operator, a, 
in row 1. Denote the remaining operators on row 1 by g'. Let 
V 2 be the stabiliser projector of the array obtained by deleting 
row 1 and column 1 from A 2 . Thus V 2 is of the form 

We then have 

V x V 2 Vx = ( ^— \ - 



' v[v' 2 v[ 



I t + sxX 

2 2 
ti 1 + SiX 1 + siX 

+ TT ^ ° ^ 



'V[g'rLK 



3. Case (c) 

Steps (c.i), (c.ii) and (c.iii) bring column 1 of A 2 in RREF 
form. Column 1 will contain two Pauli operators, o\ in row 
1, and a 2 7^ o\ in row 2. Step (civ) ensures, by suitable row 
multiplication or transposition, that a 2 is an X operator, so a\ 
is not. Let the remaining operators on rows 1 and 2 be denoted 
by g[ and g' 2 , respectively. Let V 2 be the stabiliser projector 
of the array A' 2 , obtained by deleting rows 1 and 2 and column 
1 from A 2 . Then V 2 is given by 

Thus 



V x V 2 V x 

,1 + sxX 



(- 



1 +txai ® g[ l+t 2 X®g' 2 



x(i®^)(i±il^®^) 



t+sxX 



• V'xV' 2 V'x 

1 + sxX 1 + sxX . , w ^ >/ 
+ h ^—ox ^— ® V[g'xV 2 V[ 

1+sxX t + 8l X 
1+sxX 



V'xgWx 



■tit 



112- 



axx^±^®v[g' ig ' 2 r' 2 r[ 



Sl X 



sxt 2 ^^®v[ g ' 2 r^r'x 



11 



giving 



1 1 
2~ 



V' 2 V[ 



1 1 
2~ 



7^2 a projector. This corresponds to T\ = 



with V 2 = 

1/2. This is implemented in steps (c.v) and (c.vi) through 
multiplying row 2 in A 2 by row 1 of A±, dividing T by 2, and 
subsequently deleting row 1 and column 1 in A 2 . 



4. End Game 

We still have to consider the situation where there is only 
one column left and the one where V\ is a tensor product of 
identity operators. 

The first situation is when Nl = N. In that case the 
symbols V[, g' and V 2 used in the previous subsections are 
meaningless. However, we can still make sense out of the cal- 
culations if we replace these symbols formally by the scalar 
1. Inspection of the relevant calculations then shows that at 
the very end of the algorithm, if Nl = N we have to check 
whether one of the remaining generator phases U is —1, in 
which case the states under consideration are orthogonal. That 
means both the overlap and the Uhlmann fidelity are 0, which 
we impose by setting T = 0. 

If Nl < N but V\ acts as the identity on columns Nl to N, 
V1P2P1 reduces to V 2 . This can easily be decomposed as a 
tensor product by calculating its rank R 2 (the easiest way to do 
this is by using the RREF algorithm). Thus the remaining part 
on columns Nl to N of V\P 2 V\ = V 2 is unitarily equivalent 
to 



1 +X 



II 



«)Af-_Ri-fl2 



;.=1 



5. Overlap and Uhlmann Fidelity 

In the previous subsection we have shown that V\P 2 V\ is 
equal to T times a tensor product of R\ + R 2 rank- 1 projectors 
and N — Ri — R 2 identity operators. Calculating the overlap 
and the Uhlmann fidelity is now easy. Assuming that i?i = 
K\, we have for the overlap 



F 



2 n-k! + n-k 2 Tr[7Wi] 
1 



2N-K!+N-K 2 
2~{N~K 2 +R, 2 ) rp 



(6) 



where it has to be noted that R 2 < K 2 . Since T is also a 
negative power of 2 one sees that the overlap takes values of 
either or 2~ 3 , where j is an integer between and N. 



Similarly, the Uhlmann fidelity between states pi and p 2 is 
given by 



F u = Tr[^/ A /prp2V / Pi]- 

Again we substitute the stabiliser states for their appropriately 
scaled projectors. Noting that the square root of a projector is 
that same projector gives 



VII. BIPARTITE NORMAL FORM 



(7) 



In this Section, we will modify the single-party algorithm 
CNF1 so that it can be used to reduce a stabiliser array of a 
bipartite system to a certain normal form. This algorithm will 
allow us to deduce the exact structure of this normal form, 
which is the content of Theorem^ This Theorem basically 
tells us that a bipartite mixed stabiliser state is locally equiva- 
lent to a tensor product of some number of pure EPR pairs and 
a separable mixed state. The main benefit of this normal form 
is that the entanglement of the state can immediately be read 
off from the normal form. Because of the Theorem, it turns 
out that in order to calculate the state's entanglement it is not 
necessary to actually compute the normal form completely. 
Instead, a simplified algorithm to calculate entanglement will 
be presented. 

Let us start with the statement of the normal form. 

Theorem 1 Consider a system of N qubits, separated into 
two parties, A and B, containing Na and Nb qubits, respec- 
tively. Consider a stabiliser state described by an array of K 
independent, commuting generators. 

i) By applying a suitable sequence of elementary row opera- 
tions and local elementary column ( qubit) operations, the sta- 
biliser array can be brought into the following normal form: 
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(8) 

Here, the asterisk stands for either 1 or X, and the double 
line is the separation line between the two parties. 

ii) Every pair of rows containing the XZ combinations cor- 
responds to two qubits (one from each party) being in a pure 
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maximally entangled EPR state and completely disentangled 
from the other qubits. The rows in the lower blocks of the nor- 
mal form, containing only 1 and X operators, correspond to 
the remaining qubits being in a ( in general, mixed) separable 
state. 

in) The stabiliser state described by the stabiliser array is 
locally equivalent to a tensor product of a certain number p of 
EPR pairs ^> with a separable state. For any additive entan- 
glement measure E, the entanglement of the stabiliser state is 
pE(ty). An upper bound on p is given by 



p<mm{[K/2\,N A ,N B ). 
Equality is obtained if and only if K = 2Na = 2Nb- 



(9) 



We will postpone the proof of part i) of the Theorem, the 
normal form, to the end of this Section. The proof of part iii), 
the entanglement properties of the normal form, is elementary 
and is left to the reader. The proof of part ii) is presented next. 

Proof of part ii). For convenience of notation, we first per- 
mute the qubits in such a way that the pairs of columns having 
an XZ pair in the same rows are adjacent. By eq. 0, the sta- 
biliser state corresponding to the normal form of the Theorem 



1 



2 N-K 



2p 

n 

fc=i 



fi<- 



9i 



l=2p+l 



From the specific form of the generators one sees that p can 
be written as tensor product 



P 



X 1 +z®z 



where p' corresponds to the factor n^2p+i(^ + #0/2 con- 
taining 1 and X operators only. It is a simple matter to verify 
that the factor (1 + X ® X) (1 + Z ® Z)/A is identical to the 
EPR state * = \i/j){i>\, with if; = (1, 0, 0, l) T /\/2. 

It is also simple to see that p' is a separable state. As it 
only contains 1 and X operators, it is diagonal in any basis 
where X is diagonal, and it is well-known and easy to see that 
diagonal states are separable. □ 



A. Algorithm CNFP for calculating the entanglement 

We now present an algorithm to calculate the number of 
EPR pairs in the normal form, without actually reducing the 
stabiliser array completely to that normal form. This algo- 
rithm is almost identical to algorithm CNF1, the reduction al- 
gorithm for the single-party case. 

To calculate the entanglement, the initial active region is set 
to comprise the block of elements pertaining to party A only, 
rather than the full stabiliser array, and algorithm CNFP (CNF 
for a single Party) is run on this active region. 

Algorithm CNFP is identical to CNF1, apart from the fol- 
lowing two differences: 

Step (b.vi). While in CNF1 step (b.iii) is never undone by 
step (b.v) due to commutativity of the generators, this need no 



Initial stabiliser array: 

X 
Z 

Depending on the content of column 2, do the following: 

IX 
IY 

1Z: Using a single-qubit operation, bring column 2 to &Z, then 
perform a CNOT with column 1 as target (!) and column 2 as 
control. 

XI 
Yl 

Zi : Using a single-qubit operation, bring column 2 to XI, then 
perform a CNOT with column 1 as control and column 2 as 
target. 

XX 
YY 

ZZ: Using a single-qubit operation, bring column 2 to ZZ, then 
perform a CNOT with column 1 as target (!) and column 2 as 
control. Column 2 now contains Zt. Apply another single- 
qubit operation to bring this to X 1 , and (as in the above cases) 
perform a CNOT with column 1 as control and column 2 as 
target. 



TABLE V: Required operations to eliminate all Pauli operators from 
column 2 of the stabiliser array shown above, in the various cases 
encountered. 



longer be the case here. Indeed, here we restrict attention to 
one party only, and the parts of the generators local to party A 
need not commute. Hence step (b.v) might leave Z operators 
in the leftmost column of the active region. We thus need a 
modification here: first we must check whether this has hap- 
pened and only if there are no Z operators in this column may 
Kjj and Nl be increased by 1 . Otherwise, the extent of the 
active region must stay the same. The additional Z's will then 
be treated in the next iteration of the algorithm. 

Step (civ). In step (civ), the original algorithm looked for 
an anticommuting pair in the top two rows, the presence of 
which having been guaranteed by generator commutativity. 
Here, again, this is no longer true, because the pair might be 
located in party B, which we are not allowed to touch here. 
We therefore need a second modification, to deal with the case 
that there is no such anticommuting pair. In that case, instead 
of steps (c.v), (c.vi) and (c.vii), the following operations must 
be executed. Recall that the first column has an XZ pair in 
its first two rows. This pair can now be used to eliminate all 
other Pauli operators in both the first two rows (by suitable 
single-qubit operations and CNOTs) and in the first column 
(by suitable row multiplications). TableslVland lVII contain the 
details. After that, the XZ pair can be split off from the active 
region to form part of the normal form, by increasing Kjj by 
2, and Nl by 1. 

Algorithm CNFP brings only that part of the stabiliser array 
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Initial stabiliser array: 




Depending on the content of row 3, do the following: 
X: Multiply row 1 with row 3. 

Y: Multiply row 1 with row 3, and then row 2 with row 3. 
Z: Multiply row 2 with row 3. 



subarray consisting of party B's lower right block (restricted 
to that block's columns) consists of mutually commuting gen- 
erators. By applying algorithm CNF1 to that subarray it can 
be brought in single-party normal form, consisting of X and 
1 operators only. Evidently, the row operations performed 
by the CNF1 algorithm (row permutation and multiplication) 
will also affect the corresponding rows in party A. However, 
as party A has only X and 1 operators in those rows, no Y 
or Z operators will be introduced, and the end result will also 
contain only X and 1 operators. □ 



TABLE VI: Required operations to eliminate any Pauli operator from 
row 3 of the stabiliser array shown above. 

in normal form that belongs to party A. Nevertheless, this is 
enough to read off the number of EPR pairs in the full reduc- 
tion. This will be proven below. The number of EPR pairs p is 
simply given by the number of XZ pairs in the normal form 
of party A. 

B. Proof of part i) 

By suitable modification of the Proof of algorithm CNF1, 
it can be shown that algorithm CNFP brings that part of the 
stabiliser array belonging to party A to the form as shown in 
(|8}, the columns left of the double vertical line. 

We next show that by further applying suitable column op- 
erations on the columns of party B, the complete normal form 
of (|8} can be obtained. 

Consider the first XZ pair in party A. By commutativity 
of the generators, there must at least be 1 anticommuting pair 
on the same rows in party B. By a column permutation and a 
suitable single-qubit rotation, this anticommuting pair can be 
moved to the first column of party B and be brought in XZ 
form. Using suitable CNOTs (see TablefVl the operators right 
of the XZ pair can all be brought to an 1 operator. Again 
by commutativity, the operators below the XZ pair must then 
automatically be all 1 operators. Indeed, if a row (below the 
second) contained a Pauli operator in the first column of party 
B, it would not commute with either the first row, the second, 
or both. 

One can proceed in a similar fashion with the second of 
party A's XZ pairs and party B's second column and third 
and fourth row, and so forth until all of A's XZ pairs have 
been treated in this way. 

What remains then are the rows below the horizontal line in 
(|8j. To show that the lower right block of party B in (jSJl can be 
brought to the form as advertised (i.e. containing only X and 
1 operators, as denoted by the asterisks), we note that party 
A contains no anticommuting pairs in those rows. Hence, the 



VIII. CONCLUSION 

The stabiliser formalism is a convenient tool for the study 
of entanglement properties of large quantum many -body sys- 
tems. While the stabiliser formalism provides an efficient de- 
scription of the quantum state in terms of eigenvalue equa- 
tions, it is not immediately obvious how to obtain physi- 
cal properties directly from these eigenvalue equations with- 
out explicitly having to write out the corresponding quantum 
state. In this paper we have presented, employing elementary 
tools, a number of normal forms for pure and mixed stabiliser 
states. We have furthermore provided explicit, detailed de- 
scriptions of algorithms, whose convergence we have proven, 
that allow the generation of these normal forms. Using these 
normal forms, we can compute any entanglement measure, 
overlaps between stabiliser states and various other quanti- 
ties. Detailed descriptions of the algorithms are provided that 
should make it straightforward to implement them in any pro- 
gramming language and we are able to provide MatLab suite 
of programs on request. 

These algorithms provide a firm basis for the exploration 
entanglement properties of stabiliser states and suitable gen- 
eralisations in a great variety of contexts. For example, it is 
readily seen that our approach is suitable for the efficient sim- 
ulation of systems where the initial state is a linear combi- 
nation of a polynomial number of stabiliser states. This and 
other applications will be explored in forthcoming publica- 
tions. 
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